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Integral equation formulations of electromagnetic scattering prob- 
lems have been successful for scatterers which are only a few wave- 
lengths wide. For a large scatterer, the integral equation methods 
have been uneconomical because too many basis functions are 
needed to represent the surface current on the scatterer. We use 
scattering from an infinitely long thin strip to demonstrate that the 
integral equation methods are economical for wide scatterers if the 
basis functions are properly chosen. Exact solutions of simple prob- 
lems, as well as the geometric theory of diffraction, can be used to 
suggest the appropriate functions. The results demonstrate the fea- 
sibility of the techniques used. 

I. INTRODUCTION 

Integral equation methods for solving electromagnetic scattering 
problems reduce a two-dimensional problem to a one-dimensional 
integral equation by formulating the problem in terms of the unknown 
surface currents on the scattering body or bodies. Similarly, three- 
dimensional scattering problems are reduced to two-dimensional in- 
tegral equations. 

The unknown surface currents are approximated by a sum of basis 
functions times unknown coefficients. The basis functions are usually 
chosen to be simple piecewise constant or linear functions, although 
higher-order polynomials or piecewise polynomials are sometimes 
used. The coefficients are chosen to solve (approximately) the integral 
equation. The resulting electromagnetic field exactly obeys Maxwell's 
equations and the radiation condition, but obeys the boundary condi- 
tions on the scatterer(s) only approximately. 

Even if higher-order polynomial basis functions are used, the density 
of functions on the scatterer must be at least two per wavelength. As 
a rule of thumb, five to ten functions per wavelength are typical. For 
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scatterers whose dimensions are more than a few wavelengths, the 
integral equation formulation takes too many of the usual basis func- 
tions, and has therefore been thought to be uneconomical. 

For scatterers which are many wavelengths wide, the integral equa- 
tion formulation can be economical if sufficiently good basis functions 
can be found. This is demonstrated using a two-dimensional problem, 
consisting of scattering a plane wave from a perfectly conducting, 
infinitesimally thin strip x = 0, — d/2 as y < +d/2, — oo < z < oo. This 
geometry has been the subject of considerable analysis, and approxi- 
mate solutions are available. The exact solution as an infinite series of 
Mathieu functions is also available, but is difficult to use. 

For this model problem, excellent accuracy can be obtained with 
only 17 basis functions, even for a strip more than 60 wavelengths 
wide. 

If the incident wave depends on z only as exp(i^ 2 z), the problem is 
two dimensional, and the vector electromagnetic scattering problem 
may be separated into two scalar scattering problems, one for the E- 
wave, and one for the H-wave. This is done in Section II. 

In Section III, we present integral equations to be solved for the 
currents on the scatterer. 

Section IV takes up the difficulties in solving the integral equations. 
An appropriate set of basis functions may be found by considering the 
form of the extract solution for scattering from a half-plane. (Some of 
the integrals are either singular or Cauchy principal-value integrals, 
and must be done carefully.) 

After the approximate currents on the scatterer have been found, 
the electromagnetic fields far from the scatterer can be calculated. In 
Section V, the far fields are calculated numerically and compared to 
those predicted by the geometrical theory of diffraction, which is 
accurate in the high-frequency (wavelength «: width of scatterer) 
limit. Agreement is excellent. 

II. FORMULATION OF THE SCATTERING PROBLEM 

We assume a time variation of e -WJ ' and suppress the factor through- 
out. All z dependence of fields is exp(ik z z). Then any solution of 
Maxwell's equation may be written as a superposition of two scalar 
wave equations (Ref. 1, §11.6). This may be seen as follows. Two of 
Maxwell's equations are, in rationalized mks units, 

V X E = iw/ioH, (1) 

V X H = -iweoE. (2) 

Applying the curl operator to eqs. (1) and (2), and using V-E = and 
V-H = 0, we obtain 

V 2 V+A: 2 V=0, (3) 
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where V is any of E x , E y , E z , H x , H y , or H z , and k = co V/io€o is the wave 
number. E x , E y , H x , and H y may be expressed in terms of E z and H z . 
For example, from (1), 

iupoHy = ik z E x - — — , (4) 

ox 



and from (2) 



8H 

-Uo€oE x = — - - ik z H y . (5) 



Eliminating H y and rearranging, 

8E Z dH z \ 

We let k 2 =k 2 - kl Similarly, 



E* = —2 — u\ k '-rr + a f l0 ^ry < 6 > 



i /, dE z dH 2 



By-Til*.-^?-«Po-£rl. (7) 



»-al*-S— *J- 



i / dif 2 , d# 

T2 1** -r- + W€o_ r 



H,- Ti \k,- r ± + ue - r L \. (9) 



From (3), E z obeys the scalar reduced wave equation: 

^i + £j| + A ?E , = 0, (10) 

3x dy 

as does ff*. 

The boundary conditions do not couple the E wave (E z ^0,H z = 0) 
and the H wave (jB 2 = 0, H z ¥* 0). Let n be the unit normal vector 
pointing into the scatterer at any point. The scatterer has n z = 0. Since 
E and H are zero inside a perfect conductor, the usual boundary 
conditions of continuity of tangential E and normal H are 

nx£ = 0, (11) 

n-#=0. (12) 

The first yields 

E z = (13) 

as the boundary condition for the E wave. Using (8) and (9) we see, 
since n z = 0, that 

/ dH z dE z \ ( dH z dE z \ 

n x k z — coe — — + n y \ k z — V toco — — = 0. (14) 

\ dx dy J \ dy dx ) 

INTEGRAL EQUATIONS FOR SCATTERING PROBLEMS 1895 







Since E z = for the H wave, (14) reduces to 

n.V# z = — = 0, (15) 

dn 

as the boundary condition for the H wave. If k z = 0, we obtain (15) by 
combining (6), (7), and (11). 

Since the equations and boundary conditions separate, we may solve 
two scalar scattering problems. 

Geometry of the problem 

We use standard spherical coordinates, (r, 6, <f>). Incident quantities 
are denoted by superscript i and scattered quantities by superscript s. 
The incident wave vector is 

k l = k(— cos <f> 1 sin 0\ —sin <j>' sin 9 l , —cos l ) 

= (k x , ky, k z ). (16) 

Letting r = (x, y, z), we write an incident E wave as 

—cos <£' cos 0'\ 

sin <f> 1 cos & , (17) 

sin 6' J 

(-sin <f>'\ 
cos4>' , (18) 

/ 

and an incident H wave as 

/ sin</>' \ 
Ek= Cne'H -cos <J>'j, (19) 

, /— cos <f> 1 COS ^'\ 

Hk = C// — e*- r ( -sin <£' cos ^ , (20) 

^ \ sin^' / 

where Ce and Ch are arbitrary complex constants. Any electromagnetic 
plane wave is the sum of an E wave and an H wave. 

III. INTEGRAL EQUATION FORMULATION 

To find the z components of the scattered fields, from which the 
other components may be found, we need to solve two integral equa- 
tions. The integral equation for E\ is fairly standard. First, suppose 
the scatterer is a perfect conductor with finite thickness; its cross 
section in the 2 = plane is a thin rectangle, I\ centered at the origin. 
Then for r on T, we have [Ref. 2, eq. (4.26)] 

Ei(r) = - l -\ m\k t R) *-^ ds. (21) 

4j r dn 
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Here E z = E s z + E' z , H§ ] is a Hankel function 3 representing outgoing 
waves, R is the distance from r to any point at arc length s on T, and 
n is the unit normal vector at s, pointing into the scatterer. This 
formulation ensures that the field obeys Maxwell's equations and that 
the scattered field obeys the radiation condition. 

As the thickness of the strip goes to zero, the contributions of the 
ends go to zero, and we get 

pd/2 

EUO, y, 0) = -i/A m\k t \y - y'\)J%(y') dy', (22) 

J-d/2 

where J%(y') is the sum of dE z /dn on the x > and x < sides at point 

We use E z (r) = Ce sin 0'e* kr , and define u = k t y, v = k t y', 
c = ktd/2, and y = k y /k t . Then 

C E sin 0' e* = - 4r \ H^(\u-v \)J%(v/k t ) dv. (23) 

4k ' J-c 

Finally, let Je{v) = J%{v/k t )/Ak t iCE sin 0\ and we have the integral 
equation for Je(v): 



■f 



H^(\u-v\)J E (v)dv. 



The usual integral equation for dH z /dn, eq. (4.27) of Ref. 2, is not 
useful for zero-thickness scatterers. A new integral equation, a gener- 
alization of one previously derived for thick scatterers, 4 has been 
derived by Morrison. 5 His work is directly applicable to thin bodies. 

JM y) is the difference of H z on the x > and x < sides. Letting 



jh(v) = - „.:j;v;:\^«i , w 



JUv/kt) 
4(k x /k t )CH(k/o)no)sine 

* and letting U\ and 1/2 be any values such that 

-c as Wi < u < u 2 «S c, (25) 

we have the integral equation for Jh(v): 

. rH[»(u-v) r . w A ( c Hliv-u) _ . . . 

e 'y = Jh(v) dv + Jh(v) dv 

u — v I v — u 

J-C J U-i 

"2 



+ m ] {\u-v\)J H {v)dv 



I H?{\u-v\ 

Ju, 



(26) 



+ f H\ u (\ u - v\) sgn(v - u) dJ " {v) dv 



- Jh(u 2 )HP(u 2 -u)- Jh(u{)H[ 1) (u - Mi). 
INTEGRAL EQUATIONS FOR SCATTERING PROBLEMS 1897 



Here HP is a Hankel function, 3 and sgn is the sign function: sgn(u) 
= +1 if u > 0, and -1 if u < 0. The final integral in the integral 
equation is a Cauchy principal- value integral. The equation simplifies 
considerably in form, if U\ = —c and u 2 = c, since Jh(—c) = Jh(c) = 0. 
For the numerical examples in Section VI, u\ and u 2 are taken to be 
close to u. The equation looks more complicated than if Mi = — c and 
«2 = c, but the integrals are no harder. 

Far fields 

Once the currents on the strip have been found, the field can be 
calculated at points off the strip by the usual formulas: 2 

• rd/2 

El(r) = exp(ik z z) - A m ) (k t R)J%{y') c?/, (27) 

4 J-d/2 

h C d/2 n R 

HKt) = -expffe) ^ H[ 1] (k t R) -4- <W) dy, (28) 

4 J- d/2 R 

where R is the vector from the point (x, y, 0) to the point on the strip 
at (0, y', 0). The other field components can be found from eqs. (6)- (9). 
For the far fields, R » d, the equations simplify greatly, since the 
asymptotic form of the Hankel functions can be used, 3 

H^(& = (2M) 1/ V ({ -™ /2 - ff/4) + Oif 1 ). (29) 

We let p = (x 2 + y 2 ) 1/2 . For R » d, 

R = [x 2 + (y-y') 2 f /2 = P -^+o(^j. 

Only the first two terms need be retained; the second term need be 
retained only in the exponential. Thus 

H£\k t R) = (4-Y 2 

\irktpj 

.exip[i(k tP - ktyf/p - nir/2 - tt/4)] + 0(d/p). (30) 
We use the auxiliary function F, 

F(i) = l-j e il *-~' 4 K (31) 

Finally, using the definitions of Je and Jh, changing variables to v — 
k t y\ and using a superscript /"to denote far fields, we get 



E{(r) = -exv(ik z z) C E sin 0'F(M e - ivy/p J E (v) dv, (32) 

J—c 

k x r 

H f Ar) = -exp(ik z z)C H sin O'-^-Fiktp) e - ivy/p J H (v) dv. (33) 

*' P J-c 
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The integrals may be evaluated by the methods described in Sec- 
tion IV. 

For incident waves with 0' = it/2 (k z = 0), far-field coefficients may 
be defined (Ref. 6, p. 6): 

E wave: P E = Ei(r)/C e F{k t p), 
H wave: P H - Hi(r)/C H F(k t p). 

IV. NUMERICAL SOLUTION OF THE INTEGRAL EQUATIONS 

We first give an overview of the problem, and then discuss various 
aspects of the solution method. 

4.1 Overview 

The problem has now been reduced to solving two uncoupled one- 
dimensional linear integral equations on the boundary of the scatterer. 
The usual method to solve such integral equations approximately is 
first to approximate the unknown J by a sum of basis functions, {/„} : 

J(V) - J Clnfn(v). (34) 

n-1 

The choice of the {/„} will be discussed later. Then the integrals are 
done, reducing the integral equation to 

e iyu = I aJniu), (35) 

n-l 

where I n (u) is defined by replacing J(v) by f n (v) in the integral 
equation, and performing the integrals. For the E equation, for 
example, 

In(u)=l m\\u - v\)f n {v) dv. (36) 

If the basis functions are simple piecewise polynomials, the integrals 
can be done simply by (say) the trapezoidal rule, or by Gaussian 
quadrature. For the basis functions which we will use, the integrals 
must be done much more carefully. 

For finite N, eq. (35) cannot hold for all a in (— c, c). We choose M 
points {«>}, M > N, at which to make (35) hold approximately. For 
ease of computation, the M equations in N unknowns are usually 
solved in a least-squares sense. Other approximate solutions are also 
reasonable, though more difficult to calculate. M ~ 2N is a common 
choice, and it is used for the calculations in Section V. 

4.2 Basis functions 

In scattering problems where the scatterer is not many wavelengths 
wide, a common practice is to let the basis functions {/"„} be piecewise 
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constant or piecewise linear functions, zero over all but a small part of 
the scatterer. Then the integrands of the integrals defining the {/„} 
are nonzero over only a small part of the scatterer, and evaluating the 
integral is relatively easy. Frequently the integrals are evaluated by a 
simple midpoint rule, except for the part of the integral with v near u, 
where the integrand is singular. When higher-order piecewise polyno- 
mials (splines) are used, the integrals are more difficult, but fewer 
basis functions are needed. 

We are interested in scatterers whose widths are many wavelengths, 
perhaps hundreds. We see shortly that the correct currents J have 
some components that vary as e' yu and others as e ±,u . If splines are 
used as basis functions, at least two splines per cycle of a sinusoid are 
needed, and more if reasonable accuracy is to be obtained. 

In addition, the currents near the edges of the scatterer are not well 
approximated by polynomials, and many additional basis functions 
would be needed there. 

Basis functions which closely mimic the true currents are needed, so 
that only a few functions will suffice. Of course, for most problems the 
true currents are not known. (For the thin strip, the true currents are 
known analytically, 7 but only as an infinite sum of Mathieu functions. 
The analytic solution is of less use than a purely numerical solution, 
since the individual terms in the sum are hard to evaluate, and since 
the sum converges slowly.) 

One approximate component of the true current is known, the 
physical-optics approximation to the current, or the current that would 
flow in an infinitely- wide scatterer. Far from the edges of the scatterer, 
the current should be approximately equal to the physical-optics 
current, which is proportional to e lyu . Therefore, a reasonable basis 
function is e tyu . To better approximate the physical-optics-like part of 
the current, we may also include a few (u/c)e' yu , (u/c) 2 e lyu , • • •. These 
terms are equivalent to approximating part of the current by J po (u)e tyu , 
and expanding Jpo(u) in simple smooth basis functions, since J po (u) is 
not expected to vary rapidly. 

Other approximate components of the current can be obtained from 
the exact solution to a simpler problem, reflection of a plane wave 
from a thin half-plane. Scattering from a half-plane was solved by 
Sommerfeld, for a plane wave incident at l = m/2. The solutions are 
readily accessible in Ref. 7 and in Chapter 11.5 of Ref. 1. The currents 
are conveniently written in terms of the angle a = tt/2 — <£', the angle 
between the incident wave and the positive y axis: 

J%(y) = sina<?-^ c08a 

-it.ky-n/4) 



{2 sinaG[V2A>y cos(a/2)] - iy/2/k y sin(a/2)}. (37) 



sir 
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Here G(£) is the Fresnel function, 1 

G(£) = e-* ( e*dp. (38) 



For large positive £, 



For small positive £, 



G(0 = 2|+0(r 3 )- (39) 



G(|)=i>/^' v/4 -! + 0(! 2 ). (40) 

The i/ current is 

o 

^(•y) = g-^cosa _ e'^-^Gtv^ cos(a/2)]. (41) 

•V77 

For large y, away from the edge of the reflecting half-plane 1 

J%{y) = sinae-'* cosa + e*' '■ 0[(^)" 3/2 ]. (42) 

The current is equal to the physical-optics current plus correction 
terms which oscillate and decay away from the edge. Note that the 
oscillation of the correction term is ft** rather than e** 00 * 01 ; the correc- 
tion terms correspond to waves propagating along the surface. 
Far away from the edge, the H current is 

JHy) = e- ikycoRa 

- 8eC( ^ 2) e i{hy+ " /A) {(ky)- l/2 + 0[(ky)- 3/2 ]}. (43) 

The correction to the physical-optics current decays more slowly with 
distance from the edge than does the corresponding correction to the 
E current. 

Near the edge, 

J%{y) = (2/ir) l/2 e i(ky -" /4) sm(a/2){(ky)- 1/2 

+ 4 cos 2 (a/2)(^) 1/2 + 0[(ky) 3/2 ]}, (44) 

which is infinite at the edge itself. The H current is similar, but it is 
not infinite at the edge: 

JUy) = 2( - I cos(a/2)e' ( ^-" /4) {(^) 1/2 + 0[(ky) 3/2 ]). (45) 

Suitable basis functions for the strip can be obtained by using similar 
functions at each edge of the strip. We use different basis functions 
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near the edges and away from the edges. Near the left edge, in (— c, 
—c + 8), we use basis functions 

' j. \ m-1/2 
e ilc+v) i^\ , (46) 



where m = 0, 1, 2, • ■ • for J E and m — 1, 2, • • • for </ w . The physical 
optics current is not used near the edges. Near the right edge, in (c — 
6, c), we use basis functions 

m-l/2 

(47) 



Mc-v) ' C ~ V 



Away from the edges, in (8 — c, c — 8), we use three families of basis 
functions: physical optics, decaying from the left edge, and decaying 
from the right edge. With m = 0, 1, 2, • • • , these are 



physical optics: I - ) e' 



c . 

m+l/2 



,i(c+v) 



C + V, 

m+l/2 

i(c-u) 



decaying from left edge: 



decaying from right edge: I e 

\c-v/ 

The interference between waves scattered at the two edges causes 
the coefficients to be different from the corresponding half-plane 
coefficients, but the same forms suffice. The solution is made more 
complicated by having several different kinds of basis functions, but 
only a few of each will be needed. 

The same basis functions would be obtained by taking the surface 
current calculated by geometrical diffraction theory, 7 and expanding 
the current both near and far from the edges of the strip. 

4.3 Integrals 

There are two sources of difficulty in performing any of the integrals 
to obtain any I n {u) for a given point u. [We assume none of the points 
±c, ±{c — 8) will be chosen for a u point.] First, there are singularities, 
both in the Hankel functions at u = v and in the basis functions at 
v = ±c. Second, both the Hankel functions and the basis functions are 
oscillatory, with many cycles present in some of the integrals. Many 
kinds of integrals must be done; we discuss only a few of them. 

We first consider integrals over the basis functions near the edges. 
As an example, consider 

S-C I \ 772-1/2 

m\\u-v\)e i(c ^[ C -^-\ dv, (48) 



1 
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for u > 8 - c. The Hankel function is not singular in the interval. If w 
= (c + v)/8, the integral is 



5 Hh 1) (u + c-w8)e iwS w m - l/2 dw. (49) 



Since m ^ 0, this integral is done by Gauss quadrature formulas with 
weight function w~ l/2 . A series of Gauss quadrature formulas with 
increasing numbers of points is used, until two successive values agree 
sufficiently closely. If u is within the interval, the singular part centered 
around u = v is handled separately, and the part nearer the edge is 
handled similarly to the above method. 

As another example, consider 

H$\\u-v\)f{v)dv t (50) 

H 

for some u 2 > u. (In practice we choose u 2 «s u + 1.) We let w 2 = v 
— u, obtaining 

r(u 2 -u)W 

2 wHP(w 2 )f{u + w 2 ) dw. (51) 

Jo 

The integrand is no more singular than w ]n(w). We use a complex 
version of Patterson's automatic quadrature program. 8 

The small argument behaviors of the Hankel functions are 3 

Ho l) (0 = - ln(£) + Od) + 0[£ 2 ln(fl], (52) 

77- 

HP<0--^ + m»<& (53) 

irk 

SP(Q = - ]n(& + 0(0 + OU 3 ln(|)]. (54) 

IT 

For < £ < 8, HfP(Q and BP(& are evaluated using the approxima- 
tions in Ref. 9. For example, 

BF(& = M& + iYott) = JoU) + i f |^o<e ln(jc) + ?o(&\ (55) 

Both Jo(£) and YoU) are approximated by rational functions of (£/8) 2 . 
In the example of Section V, numerator and denominator are taken to 
be quintic polynomials, giving an absolute error less than 10" in 

Jo (£) and YoU). 

For large values of £, the Hankel functions decay slowly and oscil- 
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late. 3 For £ > 8, other approximations from Ref. 9 are used. For 
example, 

BPiQ = (2M) 1/2 [P (£) + iQo(S)]e il *-* /4 \ (56) 

Both P and Q are approximated by polynomial functions of (8/£) 2 . In 
the examples of Section V, cubic polynomials are used, giving an 
absolute error less than 10 -8 . 
Integrals such as 



( 



u+8 

(1) 



m i} (v-u)f(v)dv, (57) 



are computed by the Patterson method. 

For scatterers which are many wavelengths wide, c » 1, in integrals 
such as 



£ 



m i} (v-u)f(v)dv, (58) 



the Hankel function may have many oscillations, as may f(v). If f(v) 
is a basis function decaying from the right, the e l(v ~ u) factor of the 
Hankel function just cancels the oscillations of the e l(c ~ v) factor of the 
basis function. If /( v) is a basis function decaying from the left, the 
oscillations reinforce each other. If f(v) is a physical-optics basis 
function, partial cancellation or reinforcement occurs. 

If the total oscillation in the exponential part of the integral is small, 
the integral is computed by the Patterson method. Any cancellation in 
the exponential is done explicitly. 

If the total oscillation in the exponential part of the integral is large, 
a different method must be used to be economical and to avoid too 
much error from accumulation of round-off errors. Many methods 
have been suggested for doing integrals of the type 



I 



g(v)e™dv, (59) 



where w\b — a|» 1. The most widely known is Filon's method 10 ; 
(a, b) is divided into an even number of intervals, each of length h, 
and g is evaluated at the endpoints of the intervals. The intervals are 
paired, and on each pair the three values of g are interpolated by a 
quadratic polynomial. No continuity is imposed on adjacent polyno- 
mials. The polynomials times e" 1 "' are integrated analytically. Because 
no continuity is imposed, the error is proportional to ah; for good 
accuracy u>h « 1 is needed. Thus the number of sampling points must 
increase as w | b — a \ increases. 
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Methods which do not require the number of sampling points to 
increase, involve approximating g( v) on the whole interval. A series of 
papers 11 " 16 deal with approximating g by a sum of either Legendre or 
Chebyshev polynomials orthogonal on (a, b). The polynomials times 
e 1 " are then integrated. For Legendre polynomials, the integration can 
be done analytically, but care must be taken in doing the summation 
to produce a stable method. For Chebyshev polynomials, the integra- 
tion cannot be done analytically, except in terms of infinite series, but 
no stability difficulty arises. 

Many of the particular types of g( v) that we need, decay away from 
one or both ends of the interval towards the middle, and have a 
singularity just outside the range of integration. For this type of 
function, splines are better suited than polynomials; fewer terms are 
needed for a prescribed accuracy of approximation. 

Splines have been used before for integrals like (59) 17 ; however, a 
uniform spacing of sampling points was used, which is uneconomical 
for our usual integrands. Instead of a uniform mesh, we use an 
adaptively generated nonuniform mesh. 18 The splines times e luv are 
integrated analytically. Typical values of uh, where h is a mesh spacing, 
range from 1 to 50. For c = 200 and 5 = 1, typical meshes have 20 to 
40 intervals. 

4.4 Solving for the coefficients {a„} 

By methods like those just described, each I n (u) may be calculated 
for any prescribed u. We choose M points [uj) , M > N, and require 
(35) to hold in a least-squares sense at the M points. A standard 
subroutine 19 is used to solve the equations for the {a n } . Continuity of 
the J(v) at ±(c — 8) is enforced by including extra equations and 
weighting them heavily. 20 

Apparently, no useful theory exists for choosing the { uj) . Each uj 
effectively "samples" strongly near uj, and less strongly far from Uj. It 
is clear that sampling must be done in the region where each of the 
basis functions has appreciable magnitude. Some basis functions are 
nonzero only within distance 8 from an edge, and sampling points must 
be clustered there. Other basis functions decay away from the edges, 
starting at distances larger than 8, and need sampling points not too 
far from the edges. Finally, the physical-optics basis functions are less 
localized; a good distribution of sampling points for the decaying basis 
functions should also suffice for the physical-optics basis functions. 

For the numerical examples of Section V, the sampling points in 
(— c, 8 — c) and (c — 8, c) are evenly spaced. In the left half of the 
central region, the distances of the sampling points from u = 8 — c 
increase geometrically, starting with u = 38/2 - c, and ending with u 
= 0. The points in the right half are symmetrically placed. If c is small 
enough, uniform spacing is used. 
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Table I— RMS and ERR for E and H waves, M varied 





m, 


M 




E 


wave 






H 


wave 




m r 


RMS 


ERR 


RMS 


ERR 


4 


9 


17 


1.4(- 


-3) 


5.1(- 


-4) 


3.7(- 


-4) 


3.2(- 


-4) 


5 


11 


21 


1.7(- 


-3) 


1.4(- 


-3) 


2.K- 


-3) 


7.7(- 


-4) 


7 


17 


31 


2.0(- 


-3) 


9.7 (- 


-4) 


3.6 (- 


-3) 


5.9(- 


-4) 


9 


21 


39 


2.K- 


-3) 


5.3 (- 


-4) 


4.4 (- 


-3) 


2.4 (- 


-4) 


11 


25 


47 


2.3 (- 


-3) 


1.9 (- 


-4) 


4.9 (- 


-3) 


2.8 (- 


-4) 



V. NUMERICAL RESULTS 

In this section we present selected numerical results for scattering 
from a thin, infinitely long strip. All programs were run in single 
precision on a Honeywell 6000 series computer, whose relative roundoff 
error is 1.5 x 10" 8 . 

In order to compare our results with previous work, all results are 
presented for a plane wave incident at l = m/2 and </>' = tt/4. 

We let n e by the number of edge basis functions in each edge, rid the 
number of decaying basis functions from each edge, and n p the number 
of physical-optics basis functions. The total number of unknown coef- 
ficients is N = 2n e + 2n d + n p . 

For c » 1, the geometric theory of diffraction can be used to 
compare far fields. Expansions for the E wave far field are given in 
Ref. 7, p. 200, and those for the H wave in Ref. 7, p. 218: 

Pe/c = poe + Piec~ 5/2 + 0(c- 7/2 ), (60) 

Ph/c = p 0H + Pihc~ 3/2 + 0(c~ 2 ). (61) 

We used two terms of the expansion of Ph and one term for Pe, for 
c 2b 50. The maximum value of P/c is approximately n/2/2 , at specular 
reflection, <£" = — -tt/4. We checked values of P at 33 evenly-spaced 
angles in the range —tt/4 — 4tt/c ^ <J> 8 ^ —n/4 + 47t/c. 

With c = 100, S = 1, n e = rid = 4, n p = 1. so that N = 17, we varied 
M, the number of sampling points (Table I). Letting m e be the number 
of points in each edge, and m c the number in the central region, M = 
2m e + m c . 

In Tables I to V, rms indicates the root-mean-square error in solving 





Table II — RMS and ERR for E and H waves 


, n p varied 






E wave 






H wave 




n p 


RMS 




ERR 


RMS 




ERR 


1 

2 
3 
4 


2.1 (-3) 

2.1 (-3) 
1.8 (-3) 

2.2 (-3) 




7.1 (-4) 
6.6 (-4) 
7.8 (-4) 
8.4 (-4) 


4.1 (-3) 
4.0 (-3) 
3.9 (-3) 
3.8 (-3) 




3.5 (-4) 

2.6 (-4) 
3.4 (-4) 
2.4 (-3) 
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Table III — RMS and ERR for E and H waves, 


8 varied 






E wave 






H 


wave 




s 


RMS 




ERR 


RMS 






ERR 


0.5 
1.0 
1.5 
2.0 


3.6 (-3) 
2.4 (-3) 

1.3 (-3) 

2.4 (-3) 




3.5 (-3) 
2.0 (-4) 
2.7 (-3) 
5.2 (-3) 


5.3 (-3) 
3.9 (-3) 

1.4 (-2) 
3.3 (-2) 






1.5 (-3) 
2.9 (-4) 
2.1 (-3) 
1.3 (-2) 



the over-determined linear system of equations. [Because of the two 
constraint equations which impose continuity of J at ±(c — 8), the 
equations are over-determined even for M = N.) err indicates the 
maximum absolute difference between the calculated far field and the 
geometrical theory of diffraction far field, at the 33 angles checked. 
We write 1.4 (-3) for 1.4 X 1(T 3 . The errors are not strongly dependent 
on M/N. 

With c = 100, 8 = 1, n e = n d = 4, m e = 8, and m c = 17 + 2 np (so that 
M s 2A0, we varied n p (Table II). Adding additional physical-optics 
basis functions beyond the first is essentially useless. 

Table IV — RMS and ERR for E and H waves, n B varied 







E 


wave 






H 


wave 




n e = rid 


RMS 






ERR 


RMS 






ERR 


2 

3 
4 
5 


2.9 (-2) 

5.6 (-3) 
2.4 (-3) 

1.7 (-3) 






5.3 (-3) 
1.0 (-3) 
2.0 (-4) 
1.8 (-4) 


1.4 (-1) 
2.9 (-2) 
3.9 (-3) 
7.6 (-4) 






2.1 (-2) 
1.3 (-3) 
2.9 (-4) 

2.2 (-4) 



With c = 100, n e = n d = 4, n p = 1, m e = 8, and m c = 25, we varied 8 
(Table III). The best value is approximately 1. Similar results are 
obtained for c = 50 and c = 200. 

With c = 100, 6=1, and n p = 1, we varied n e and n d , keeping m e = 
2n e and m c = 4m d + 3, so that M ^ 2N (Table IV). This illustrates 
convergence with increasing N. 

Finally, with 8 = I, n e = n d = 4, n p = 1, m e = 8, and m c = 19, we 
varied c (Table V). Convergence is similar for c = 50, 100, and 200. 





Table V — RMS and ERR for E and H waves, c varied 




E wave H wave 


c 


RMS ERR RMS ERR 


50 
100 
200 


1.7 (-3) 3.0 (-4) 3.6 (-3) 7.0 (-4) 
2.4 (-3) 2.0 (-4) 3.9 (-3) 2.9 (-4) 

1.8 (-3) 0.3 (-4) 3.5 (-3) 2.2 (-4) 
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